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We study the finite-temperature behavior of the Lipkin-Meshkov-Glick model with a focus on 
correlation properties as measured by the mutual information. The latter, which quantifies the 
amount of both classical and quantum correlations, is computed exactly in the two limiting cases 
of vanishing magnetic field and vanishing temperature. For all other situations , numerical results 
provide evidence of a finite mutual information at all temperatures except at criticality. There, it 
diverges as the logarithm of the system size, with a prefactor that can take only two values, depending 
on whether the critical temperature vanishes or not. Our work provides a simple example in which 
the mutual information appears as a powerful tool to detect finite-temperature phase transitions, 
contrary to entanglement measures such as the concurrence. 
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I. INTRODUCTION 

The theory of entanglement, originally developed in 
the field of quantum information theory, has recently be- 
come a very valuable tool to describe the phase diagram 
and the properties of strongly-correlated quantum many- 
body systems' . For example, the study of the scaling of 
the entanglement entropy as a function of the system 
size has become the standard tool to detect the central 
charge in critical quantum spin chains"' ', and the entan- 
glement spectrum in topological insulators is providing 
a very nice characterization of the spectrum of the cor- 
responding edge modes'''. Also, second-order quantum 
phase transitions are characterized by the emergence of 
a large amount of entanglement. 

Almost all of those results were obtained in the zero- 
temperature regime, i.e., for ground states. In that 
regime, all correlations are quantum correlations, and 
the subtleties of mixed-state quantum entanglement are 
avoided. From the point of view of quantum systems at 
finite temperature, interesting questions about separabil- 
ity, entanglement cost and entanglement of distillation 
could in principle be posed, but those concepts do not 
seem to have any clear operational meaning in the con- 
text of a quantum many-body system as there is no clear 
separation of the different degrees of freedom in that case. 
Indeed, at finite temperature, it seems to make much 
more sense to look at the total correlations in the system, 
without trying to distinguish classical from quantum de- 
grees of freedom: it is precisely the fact that (quasi-) long 
range correlations emerge in the many-body system that 
makes them so interesting, and whether those correla- 
tions are quantum or classical is a question that would 
anyway depend on the way the system is partitioned and 
on the definition of entanglement. 

The most appealing quantity that measures the to- 
tal amount of correlations in a quantum system is the 
mutual information. For pure (zero-temperature) states 
it reduces to twice the entanglement entropy, and effec- 



tively quantifies the amount of (Shannon) information we 
acquire about the configuration of part of the system by 
measuring another part. For instance, in the case of an 
Ising ferromagnet at zero temperature, the mutual infor- 
mation between two parts of the system is one bit, as the 
measurement of the spins in one part reveals one bit of 
information (up or down) about the configuration of the 
spins in the other part. In such purely classical cases, the 
behavior of the mutual information as a function of the 
temperature and the system size is very intriguing and 
has been discussed in Refs. (i-lO. The really nice features 
of the mutual information are that: (i) it provides a nat- 
ural way of extending the study of entanglement entropy 
to the case of finite temperature; (ii) it has a nice oper- 
ational interpretation; (in) it provides a measure of all 
correlations in a quantum many-body system, not just 
the two- or three-body correlations. This last point is of 
special importance when studying collective models such 
as the Lipkin-Meshkov-Glick (LMG) model''"'', where 
there is no clear notion of locality. 

The goal of this paper is to initiate the study of 
the behavior of the mutual information in quantum 
many-body systems at finite temperature following re- 
cent studies""''. Questions that will be addressed in- 
clude, e.g., the study of how the scaling of the mutual in- 
formation depends on the temperature, and how precisely 
the mutual information allows one to identify phase tran- 
sitions. We will focus our attention on the LMG model, 
because it is a quantum many-body system that exhibits 
a nontrivial phase transition with true many-body cor- 
relations, while still being simple enough such that the 
calculation of the mutual information is tractable. 

We will show that the mutual information yields pre- 
cise information about the location of the phase transi- 
tion. An intriguing difference in its finite-size scaling is 
obtained between the zero versus non-zero temperature 
case. 
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II. MODEL 

A. Hamiltonian and symmetries 

Let us consider a system of N fully-connected spins 
1/2 whose Hamiltonian reads 



H 



' N 



hS, 



(1) 



where Sa — \ 'Ylij=i'^'j the total-spin operator along 
the a — X, y, z direction, 



cr" being the usual Pauli ma- 
trix at site j. Without loss of generality, we suppose 
that h ^ 0. The Hamiltonian (1) can be seen as a ferro- 
magnetic transverse-field Ising model with infinite-range 
interactions (see Ref. 18 for a recent review). It is also a 
special case of the LMG model^^"^'* whose ground-state 
entanglement properties have triggered much attention 
during the last years^''"-''. 

As can be seen easily, this Hamiltonian preserves the 
total spin, i.e., [H,S'^] — 0. This conservation rule, to- 
gether with the so-called spin- flip symmetry [H, Y[j (^j] ~ 
0, allows one to investigate this model both analytically 
and numerically. Recently, the exact spectrum has been 
obtained in the thermodynamical limit, in the maximum 
spin sector s = N/2 where the ground state lies"''""". Un- 
fortunately, a complete derivation of the full spectrum (in 
arbitrary s sectors) that is crucial for the understanding 
of the finite-temperature properties is still missing. 



B. Phase diagram 

In order to get a first feeling about the model, let us 
start with a discussion of the finite-temperature phase 
diagram. It can be obtained using a standard mean- 
field approach (see for instance Refs. and : ) that we 
briefly describe here. We introduce the order param- 
eter (proportional to the average magnetization along 
the x-direction) nix — (cj) = jj{Sx), and write that 
a J = rrix + (crj — rrix). Then one can plug this ex- 
pression into the Hamiltonian, neglect fluctuations in- 
volving (crj — rux)^ and keep only dominant terms in a 
1 /iV-expansion. This yields an effective Hamiltonian of a 
large spin S = {Sx, Sy,Sz) in an effective magnetic field 
^cff = ("^x,0, /i), that reads 



^^cfi' — —hcH ■ S, 



(2) 



up to an unimportant constant shift of N^^. One still 
has to ensure that nix coincides with the thermodynam- 
{Sx). Therefore, one needs the partition 



ical average 
function 



N 



Zcs = Tr (e" 



2 cosh 



/3\/ 



mt 



-, N 



(3) 
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FIG. 1. Density plot of the order parameter m in the {h,T) 
plane for a finite-size system (A'^ — 256). The black line in- 
dicates the transition line Tc{h) whose expression is given in 
Eq. (6). 



where we have used the inverse temperature /? 
Noticing that 



nix = Tr 



2Sx 



2 OlnZ, 



cff 



one obtains that either tUx = or nix satisfies 



y^^^ITT^ = tanh ( ^4±^ 



2T 



l/T. 
(4) 

(5) 



The above equation has a solution only for T < Tc{h), 
where the critical temperature reads 



Tc{h) = 



2tanh"^(/i) 



(6) 



This critical temperature, which is always smaller than 
1/2, is represented as a black line in Fig. 1. We thus find, 
as expected, a disordered (paramagnetic) phase at high 
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FIG. 2. Order parameter m for ft = 1/2 as a function of tem- 
perature for different values of A^. The critical exponent is 
/3mf = 1/2. The black line indicates the transition tempera- 
ture rc(l/2) ~ 0.46. 
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Along the same line, one can compute the heat capa- 
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FIG. 3. Density plot of the heat capacity per spin c in the 
{h, T) plane for a finite-size system {N = 256). The black line 
indicates the transition line Tdh). 



temperature T > (h) with a vanishing order parameter 
nia: and an ordered (ferromagnetic) phase for T < Tc{h) 
where rrix acquires a finite value. One can also easily 
show that, in the vicinity of Tc{h), the order parameter 
vanishes as nix ^ (^c — T)^"'' with the standard mean- 
field exponent /3mf = 1/2 (see for instance Ref. 30). 

The validity of the mean-field approach can be 
checked numerically, provided one considers the alter- 
native quantity m = jf\/JS^. Indeed, there is no 
symmetry breaking for any finite value of N, so that nix 
always vanishes in numerical simulations. As the system 
is fully connected, one furthermore expects that and 
TO^ become equal in the thermodynamical limit. Exact 
diagonalization results are compared with the mean-field 
prediction in Figs. 1 and 2 where an excellent agreement 
can be observed. 
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FIG. 4. Heat capacity per spin c for ft = 1/2 as a function 
of temperature for different values of TV. The discontinuity 
at T'c(l/2) ~ 0.46 (black vertical line) is associated with the 
critical mean- field exponent Qmf = 0. 



city defined as C = where U is the average energy 



U = — Tr i/efte 



d In Zes 
d(3 ■ 



(7) 



The heat capacity per spin thus reads 



iV 



InZoff 
N 5/32 ' 



(8) 



which can be straightforwardly evaluated using Eq. (3). 
As can be seen in Figs. 3 and 4, analytical results stem- 
ming from the mean-field approach are still perfectly con- 
sistent with the numerical data. Note that one recovers 
the (well-known) fact that the heat capacity displays a 
jump at the transition, so that C \T - Td"""'' with 
"MF = 0. 



C. Finite-size corrections 

Mean-field results obtained in Sec. II B are valid in the 
thermodynamical limit {N — > oo). However, it is inter- 
esting to pay attention to finite-size corrections that con- 
tain valuable informations. These can be studied quan- 
titatively since, contrary to usual quantum spin systems, 
the symmetries of the LMG model discussed in Sec. II A 
allow one to perform numerical computations for rela- 
tively large system sizes (at least N ~ 10'^ for all quanti- 
ties considered in this work). 

At zero temperature, these finite-size corrections 
have been shown to be very sensitive to the 



transition 



20,21,30,31 



. For instance, corrections to the ther- 



modynamical value of the order parameter m behave as 
iV~^/^ at the critical point and as N~^^'^ otherwise. At 
finite temperature, we found numerical evidence for a 
similar behavior for to away from criticality [more pre- 
cisely, TO^ = nix + a/N + 0{1/N^)] but, on the critical 
line Tc{h) > 0, to rather seems to vanish as iV"^/^ as 
can be seen in Fig. 5. Unfortunately, methods developed 
to analyze the zero-temperature problem ' cannot be 
used at finite temperature. This is mainly due to the fact 
that, for T ^ 0, one needs to consider all spin s sectors 
and not only the maximum s sector where the ground 
state lies. However, for h ~ (no quantum fluctuations), 
one can compute the leading finite-size corrections of all 
quantities at and away from criticality (Tc(0) = 1/2). In 
this special shown in Sec. IIIB, to is found to 

vanish as N~^^'^ giving a first indication that the non- 
trivial finite-size behavior at Tc{h) > might be inferred 
from the zero-field problem. Although one can argue 
that quantum fluctuations are irrelevant for the study of 
critical behavior at finite-temperature, it is striking to 
recover these behaviors from the trivial h = case. 
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FIG. 5. Logj-logj plot of m vs N for h = 0, 1/2, and 1 
computed at the critical temperature Tc{h) and at high tem- 
perature T — 0.7 where nix — 0. Red lines are straight lines 
with slopes —1/2, —1/3, and —1/4 (see the text). 



III. MUTUAL INFORMATION 
A. Generalities 

We now turn to the main focus of the paper, namely 
the analysis of the finite-temperature correlations. As 
already mentioned, ground-state entanglement proper- 
ties relevant for the zero-temperature problem have been 
studied extensively Note also that several quanti- 
ties based on fidelity have also been proposed to capture 
the zero-temperature transition in the LMG model ' 
but they do not measure entanglement properties. By 
contrast, there have only been a few investigations 
concerning entanglement in the LMG model at finite 
temperature^"^' * ' mainly because entanglement measures 
for mixed states are rare and hardly computable. 

This last decade, many studies have shown that en- 
tanglement measures could be used to detect quantum 
(zero-temperature) phase transitions (see for instance 
Ref. j for a review). Here, following recent works in 
two-dimensional quantum spin systems^^^^'', we wish to 
investigate the same problem at finite temperature. Con- 
cerning the LMG model, the concurrence ' ' (which char- 
acterizes the entanglement between two spins in an arbi- 
trary state) is up to now the only entanglement measure 
that has been computed both at zero^"""^ and at finite 
temperature . However, as can be seen in Fig. 6, this 
quantity is quite surprisingly insensitive to the transition 
line at finite temperature. 

Alternatively, one may consider the negativity '^ that 
is also a good measure for mixed states and that is not re- 
stricted to two-spin correlations. This measure has been 
successfully used to study the quantum phase transition 
in the LMG model'^'' and other fully-connected models . 
It is however already difficult to compute at zero tem- 
perature, so one can expect troubles to obtain it at fi- 
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FIG. 6. Density plot of the concurrence in the {h, T) plane 
for a finite-size system (A'^ — 256). The black line indicates 
the transition line Tc{h). 



nite temperature for systems of decent sizes (analytical 
expressions in the thermodynamical limit being out of 
reach). To our knowledge, apart from these two mea- 
sures (concurrence and negativity), there exists only one 
reliable quantity susceptible to be used for mixed states, 
namely the mutual information, whose definition is based 
on entanglement entropy. 

To obtain informations about many-body correlations 
of a system, the entanglement entropy is usually the mea- 
sure of choice for pure states. It is defined as the entropy 
of the reduced density matrix pA obtained, after parti- 
tioning the system into two parts A and B, by tracing 
out the density matrix p over subsystem B : pA = Tt^bP- 
For pure states, the same entropy is obtained if one ex- 
changes A and B, although the result obviously depends 
on the bipartition. For mixed states, however, the entan- 
glement entropy is dominated by a contribution related 
to the total entropy of the system under consideration 
(which vanishes for pure states). Therefore, the natural 
thing to do is to subtract this total entropy in a suitable 
way. This leads to the following definition of the mutual 
information 



I{A,B)^£a+£b-£ai 



(9) 



where the von Neumann entropy £a = — Tr(p^ log2 pa) is 
the quantum generalization of the Shannon entropy. The 
Shannon entropy provides an asymptotic description of 
the properties of a system (a probability distribution of 
states) , in the sense that it describes the average number 
of bits needed to encode a state of the system. This aver- 
age is taken over all the states of the system, using their 
respective probabilities of occurrence. In other words. 
Shannon entropy corresponds to the average amount of 
randomness that can be extracted from the system. Note 
that it is also possible to work with Renyi entropies de- 
fined as £a{k) = jz^ log2 Trp^^", but let us simply un- 
derline that when k — >■ 1, they become identical to the 
aforementioned Shannon/von Neumann version that is 
considered in the following. 
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FIG. 7. Density plot of the mutual information T in the {h, T) 
plane for a finite-size system (A'^ — 256). The black line indi- 
cates the transition line Tc(h). 



Before discussing the results, it is important to under- 
stand that mutual information is not a measure for the 
"quantumness" of a system: it captures both classical 
and quantum correlations alike. To make the distinc- 
tion between both, one could rather study the so-called 
quantum discord''^'"'". However, such a measure is clearly 
artificial as the quantum discord of two states arbitrarily 
close to each other can vanish for one of them and be 
positive for the other one"'\ Such mixed-state entangle- 
ment measures were constructed from the point of view 
of delocalized systems, and do not necessarily make sense 
in the context of condensed-matter systems. 

In the present work, we shall not try to distinguish be- 
tween quantum and classical correlations and will stick 
to computing the mutual information. Indeed, such a 
quantity should be sensitive to strong fluctuations of 
the many-body correlations arising at a phase transition. 
Numerical results are displayed in Fig. 7. Following a ray 
in the (/i, T) plane of this figure, one gets that I{A, B) 
is almost equal to 1 near the origin, grows as one goes 
away from the origin, reaches a maximum before transi- 
tion line, and then decreases to zero. Assuming that the 
maximum is located at the transition line in the thermo- 
dynamical limit (see Sec. Ill C for a confirmation of this), 
one might thus argue that I{A, B) detects the phase tran- 
sition at least as well as usual thermodynamic quantities 
such as the order parameter or the heat capacity dis- 
cussed in Sec. II B. However, one should keep in mind 
that computing the mutual information is not as easy as 
computing an order parameter, even in the simplest case 
{h — 0) that we shall now discuss. 



B. Exact results for /i = 

Let us consider the classical (quantum fluctuation free) 
case h — Q for which the Hamiltonian (1) simplifles to 



H 



' N' 



(10) 



Its eigenstates are those of Sx and can be chosen as sep- 
arable states \p, i) , with p spins pointing in the —x di- 
rection and {N — p) spins pointing in the -\-x direction 

, N 

The variable i allows to distinguish between the 



P 

such states which are degenerate and have eigenenergy 

To obtain the mutual information, one flrst needs to 
compute the total entropy of the system at finite tem- 
perature and, hence, the partition function 



N 



Z(/3) = Tre-^^=^ 

p=0 



-PEp 



(11) 



where as before, /? — 1/T denotes the inverse tempera- 
ture. Despite the apparent simplicity of H, it is difficult 
to compute Z analytically for arbitrary finite N . Thus, 
in the following, we shall focus on the physically relevant 
large- A'^ limit to analyze the thermodynamical limit and 
its neighborhood. Computations get easier in this limit 
because the discrete sum in (11) can be approximated by 
an integral with, as is well-known, an error of the order 
l/N"^. This is fortunate since, in order to compute the 
mutual entropy in the present problem, we will need to 
compute the first correction to the infinite N limit. 

The very first step in the calculation is to perform a 
large- expansion of the binomials 



TV! 



r(A^ + i) 



p\{N~p)\ T{p + l)T{N -p+iy 



(12) 



Using the series expansion of the Euler-F function at 
large N (or the Stirling formula) and skipping terms of 
relative order l/N^, one gets 



1 



ttN Vl - 4e2 

^-N[{l-e) log((i-£) + (i+£) log(i+£)] ^ 



1 - 



1 S + ie^ 
iV12(l-4e2) 



+ 0{1/N^) 



(13) 



where we set e = jj — -^^ Then, replacing the sum by an 
integral in Eq. (11) gives 



with 



->N 



2 ds 

Vl - 4£2 
1 3 + 4£2 



iV12(l-4£2) 



£ 1 log I - - £ 



0{1/N^) 



,(14) 



(15) 
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Furthermore, the exponential term e"^''''^'^^ allows us to 
extend the integration range to M. Indeed, the effect 
of this extension consists in exponentially small terms 



(cx 



) to which we anyway do not have access. 



The resulting integrals can be evaluated using the stan- 
dard Laplace's method (also known as saddle-point or 
stationary-phase approximation) though one has to take 
care of computing the subleading corrections . 

At and above the critical temperature (/3 ^ /3c (0) = 2), 
the minimum of f{e) is found for e = and the result of 
the saddle-point approximation reads 



< 2) = 2 

and 

Z(/3 = 2) = 2 



N 



2-/3 



1 fl2 

1 - — h Od/N'^) 

4(2-/3)2 -^^yy ) 



(16) 



jv_i3i/47Vi/4r(l/4) 



r(3/4) 



2807V 20^/3iV3/2^(l/4) 



^ ^ 2V3r(3/4) 

5y]vr(i/4) 

0{1/N'] 



(17) 



Note that, in this approach, Z{f3 = 2) does not coin- 
cide with lim^_j.2 Z(J3 < 2) (which is divergent) since the 
integrals stemming from the stationary-phase approxi- 
mation is not gaussian anymore at criticality. Further- 
more, we see that Z{f3 — 2)/2^ diverges with a non- 
trivial finite-size scaling exponent, as 7V^/^. This result, 
which is obtained by directly computing Z((5 — 2), is 



consistent with Eq. (16). Indeed, adapting the finite-size 
scaling argument developed in Refs. 2U and 21 for the 
zero-temperature problem, one can write (for /3 smaller 
than 2 but clos e to 2) Z( /3 < 2)/2^ = z^f[N{2 - (3)^] 
where Zoc — \/2/(/3 — 2) is the thermodynamical limit 
value of Z{/3 < 2)/2^, and / is a scaling function. As 
this quantity cannot be singular for finite N, one must 
have f{x) ^ x^l^ so that the singularity at /? = 2 disap- 
pears : Z(/3 < 2)/2^ - (/3-2)-i/2[7v(2-/3)2]i/4 _ j^ij^^ 

In the low-temperature phase (/3 > 2), •^[e] has two 
symmetric minima the position of which can only be 
computed numerically. Once determined, one can still 
perform the Gaussian integral resulting from the second- 
order Taylor expansion around these minima, which 
yields a "numerically exact" result in the infinite- A^ limit. 
As a consequence, we only give analytical expressions for 
/3 ^ 2 but we also show the exact numerical results for 
/? > 2 in the figures. 

To compute the entropy in the thermodynamical limit, 
one needs to compute the internal energy of the system 
using the same approximation (13) of the binomials. This 
quantity is defined as 



C/(/3) = Tr(pH) 



91nZ(/3) 



(18) 



where p = is the thermal density matrix. For 

/3 < 2, it can be obtained from Eq. (16) but, as previously, 
special care must be taken to deal with the critical case 
for which one can show that 



C/(/3 = 2) = - 



\/3]vr(3/4) 

2r(i/4) 
ior(i/4)4 + 



1 - 



2V3r(3/4) ^ 

5Vivr(i/4) 

32r(l/4)2r(3/4)2 



12 



r(i/4)2 + 7r(3/4)^ 

175A^r(l/4)2 
504r(3/4)4 



875^Ar3/2r(l/4)3r(3/4) 



0(1/A^2 



(19) 



so that the internal energy diverges as N^l"^ . 

Given that, t/(/3) = (U) = -^(5'2), this latter result 
directly implies that the order parameter m = (S'2) 

vanishes at criticality as N^^l^, as was already noticed 
in Fig. 5. 

The entropy can finally be computed (analytically for 
/3 ^ 2 and numerically for /3 > 2) since 

ZAB{^i) = -Tr(plog2 p) - log2 Z(/3) + (20) 

log 2 

Computing the mutual information is another game, 
the most difficult part of which lies in the derivation of 
the large- A^ behavior of the partial entropy 

£A(/3) = -Tr(pAlog2p^), (21) 

where pA = Tr^/O is the reduced density matrix. It is 
thus mandatory to find an expression for pA- With this 



aim in mind, let us split the system into two parts A 
and B containing Na and Nb = N — Na spins respec- 
tively. Next, let us decompose the eigenstates of H as 
\p,i) = \pA,iA)A ® \PB,iB)B withpA+PB =_pand where 
\pj,ij)j denotes a state of subsystem j = A,B that has 
Pj spins pointing in the —x direction. Variable ij can 

take ^ ^ values and depend on which spins point in 

the —X direction. This decomposition allows one to write 
the reduced density matrix as 
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PA = ^Trs^e 



1 



JV ^ N : 



E 



(22) 

(23) 

b(pb, isl, 

(24) 



We introduced the quantity 

1 .r- / iVR \ ^JV/3(^ + ^ 



Nb 
Pb=0 



Pi3 



RiPA) 

The partial entropy then simply reads 

R{pA)log2R{pA)- 



(25) 



Na 

PA 



(26) 



Following the same line as for the partition function cal- 
culation and replacing the binomials by the same form 
as (13), one can still use Laplace's method to obtain the 
large-iV behavior of the partial entropy. Nevertheless, 
things are a bit more involved since one now has to deal 
with a double sum and thus with a two- variable integral. 
After some algebra, one gets 



£aW <2)^tN 



1 



+ \og 



2 log 2 [ 2 
■ 2-/3 



(27) 



and 



SaW = 2) = tN- 



2-/3(l-r) 



73 r(3/4) 1 
log2r(l/4) ^4 



0(l/iV), 



log2 N + 0{N° 



(28) 

where we have set Na/N = r, and Nb/N = 1 — t. 

Finally, noting that £b is obtained from £a by ex- 
changing TO (1 — r) and using (20), one obtains the 
following expressions of the mutual information 



I(/3 < 2) 



and 



1 f [2-l3r][2-(3(l-r)] 

2 [ 2(2 - 13) 



I{I3 



■log2iV + 0(7V°). 



Oil/N), 
(29) 

(30) 



These expressions together with the exact numerical 
results for /3 > 2 (in the thermodynamical limit) are 
compared to finite-size numerical results in Fig. 8 (left 
column). In Fig. 8(a), we display the mutual entropy 
as a function of T for several values of N and r = 1/2. 
The excellent agreement with Eq. (29) can be observed 
at large N. In Fig. 8(d) we provide a numerical check 



of Eq. (30) that predicts a logarithmic divergence of the 
mutual information at criticality, with a nontrivial pref- 
actor 1/4. Note that in the present case {h = 0), it is 
possible to perform a numerical study for large system 
sizes (up to = 2^^ here) that allows one to observe 
how this divergence arises. 



C. Numerical results for h > 

Let us now turn to the case h for which the spec- 
trum, and hence I can only be computed numerically. 
Details of the algorithm we developed to compute I can 
be found in Appendix A. In Figs. 8(b) and (c), we dis- 
play the behavior of the mutual information as a function 
of temperature for h — 1/2 and h ~ 0.9999 (such that 
Tc = 0.1), and for r — 1/2. For h — 1/2, and as already 
observed for ft, = in Fig. 8(a), I reaches a finite value 
away from criticality whereas it increases with the system 
size at T^. Although this is less obvious for h ~ 0.9999, it 
is however still true. Indeed, the apparent increase of I 
away from the critical temperature observed in Fig. 8(c) 
is a finite-size artifact as can be inferred from the zero- 
temperature case. For T = and ^ h < 1, one has 
I = 2£a — 1, since £b = £a (for all r) and £ab = 1 in the 
broken phase (Sab ~ in the symmetric phase h > 1). 
Then, using the expression of £a computed in Refs. 23 
and , one gets, in thermodynamical limit, I ~ 6.1113 
for h ~ 0.9999 whereas for N = 512 we find 1 ~ 1.9378. 
It is thus clear that, in this case, the asymptotic value is 
still far from being reached for A^ = 512. By contrast, 
for h = 1/2, one gets I ~ 1.0286 in the thermodynamic 
limit and I ~ 1.0291 for N = 512. 

To investigate the divergence of the mutual informa- 
tion at criticality in more details, we computed 1{Tc) for 
increasing system sizes. Results shown in Figs. 8(d,e,f) 
confirm that I{Tc) diverges logarithmically with N . Fur- 
thermore, it seems that the behavior of I for < /i < 1 is 
the same as for /i = 0. This can be seen in Fig. 9 where we 
have superimposed data shown in Figs. 8(d,e,f). Thus, 
although numerically reachable sizes for < ft < 1 (up 
to A^ = 2^^ here) are not as large as for ft = 0, we are 
led to conclude that X{T^) - \ log^ A^ for all < ft < 1. 
Actually, if one reminds that a similar result was found 
for the order parameter m (see Fig. 5), it is tempting 
to conjecture that, for all physical quantities, the lead- 
ing finite-size corrections at finite temperatures can be 
infered from the zero-field problem. In other words, the 
system behaves classically for all nonvanishing tempera- 
tures. 

However, as already noted for the order parameter 
(see Fig. 5), things are dramatically different at zero 
temperature. Indeed, in the purely quantum case (no 
thermal fluctuations), a phase transition occurs at 
ft = 1 and the mutual information can be computed 
straightforwardly, as discussed above. 

Using the fact that £a ^ h log2 A^ for ft = 1 ' , and 
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FIG. 8. Upper row : mutual information I as a function of T for r = 1/2, various sizes, and for h = (a), h — 1/2 (b) 
and h ~ 0.9999 (c). The black vertical lines indicate the critical temperatures. Lower row : scaling of the critical mutual 
information (evaluated at T = Tc(/i)) as a function of N, for r = 1/2, and for ft = (d), h = 1/2 (e) and h ~ 0.9999 (f). The 
insets show the behaviour of the slope between two consecutive data points as a function of l/log2 A*'. Red lines indicate the 
large- behavior I ~ | logj A^. 



that Eab = at the critical point (the ground state is 
unique), one gets 

I(/i = l,T = 0)~ ilog2 7V, (31) 

in contrast vi^ith the j log2 N behavior found for 
0<h<l. 




2 4 6 8 10 12 14 16 
FIG. 9. Superposition of the curves of Figs. 8(d,e,f). 



IV. CONCLUSION 

We have studied the finite-temperature LMG model 
by focusing on the mutual information I that measures 



both classical and quantum many-body correlations in 
the system. As already pointed out in other recent 
studies''''^''"^'^ we believe that mutual information is a 
good detector of finite-temperature phase transitions. 
This quantity could turn out to be especially useful in 
systems where no obvious or simple order parameter can 
be found (for instance in topological phase transitions). 
However, this quantity may reveal difficult to compute 
even in a simple fully-connected model such as the one 
considered here. 

In the present work, we managed to get analytical re- 
sults in the simplest (purely classical) situation only but 
we obtained numerical results in the whole parameter 
range for relatively large system sizes (up to = 2^*^ 
spins). As already observed for the (zero-temperature) 
entanglement entropy^''"', we found that I is always 
finite except on the transition line Tc{h) where it di- 
verges logarithmically with the system size A'^. A com- 
plete study of the purely classical model allowed us to 
conjecture that for ^ /i < 1 it diverges as j log2 A^ 
whereas for /i = 1, 1 behaves as | logj N at criticality. 

We hope that our results will motivate further stud- 
ies of the mutual information (and more generally on 
many-body correlations) in related models since many 
questions remain open. First of all, let us underline 
that the model studied here is only a special case of 
the LMG model where spin-spin interactions are strongly 
anisotropic (Ising-like) . Thus, investigating the influence 
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of an anisotropy parameter would be valuable especially 
since for isotropic interactions (XX-like), the full spec- 
trum is known exactly (see for instance^^) and, as already 
seen at zero temperature"', one might expect a different 
behavior of I at the critical point. In addition, when 
competing (XY-like) interactions are present, there ex- 
ists a special point in the parameter space where the 
ground-state is a separable state, the so-called Kurmann 
point^ '. It would be definitely instructive to analyze 
many-body correlations near this point at finite temper- 
atures. 

Second, natural extensions of fully-connected models 
have been recently analyzed at zero temperature, reveal- 
ing a rich variety of quantum phase transitions '" among 
which some are of first order. The framework detailed in 
App. A clearly allows for a direct computation of the mu- 
tual information in these models"'^'. From that respect, it 
would also be interesting to study the finite-temperature 
Dicke model that shares many entanglement features 
with the LMG model, as discussed in Ref. 24. However, 
the presence of a bosonic mode coupled to a set of two- 
level systems may give rise to interesting challenges if one 
aims at computing the mutual information between the 
field and the atoms (in particular the procedure described 
in App. A will not be sufficient). 



Appendix A: General algorithm to compute the 
mutual information in spin-conserving Hamiltonian 

The Hilbert space of a system consisting of N spins 
1/2 has dimension 2^ which usually limits the study of 
such systems to a number of spins of a few tens. This 
restriction can be overcome in the LMG model, because 
its Hamiltonian commutes with the total spin operator: 
[H, S'^] = 0. This symmetry ensures that the Hamilto- 
nian is block diagonal and that one can study each block 
of fixed s separately. The dimension of the latter is 2s + 1 
(hence at most linear in N) which allows one to study 
larger systems (A'^ ~ 10"^ was reached in this paper). 

As is known from the theory of addition of angular 
momenta, there are distinct ways of obtaining a total 
spin s when combining N spins 1/2, where reads 



N/2 - s) yN/2 - s 

2s + 1 / N 
N/2 + S + 1 \N/2 - s 



(Al) 



As a consequence, each energy level in the spin s-sector is 
(at least) d^ times degenerate. Of course, one can check 
that 



S 

E 

s = [S] 



'2s + l)d 



(A2) 



where S = A^/2 is the maximum spin, and [S] is the 
minimum spin, namely [S] = if A^ is even and [S] — 1/2 



if A'' is odd. The fact that the sum over s runs over 
integers (half-odds) if A^ is even (odd) is implicit. 

The goal of this appendix is to explain how to use the 
spin-conserving property to compute the reduced density 
matrix from which the mutual information is extracted. 
Additional symmetries, such as the spin-flip symmetry 
[H, cr|] = present in the LMG model, will not be 
considered thereafter in order to give general results for 
any Hamiltonian satisfying [H, S'^] = 0. Of course, ad- 
ditional symmetries can be implemented to fasten the 
algorithm. 

From the discussion given above, such a spin- 
conserving Hamiltonian can be block-diagonalized and 
written as 



H 



(A3) 



s=[S] i=l 
S s 

= E E E E ht\^\s,m% .{s,ml (A4) 

s=[S] i=l m= — s 7n' = — B 

where i labels the d^ degenerate subspaces of spin s, and 
where the notations of the sums over s have already been 
introduced above. Furthermore, we introduced eigen- 
states \s,m)i of operators and with eigenvalues 
s(s -|- 1) and m respectively. It is worth noting that ma- 

trix elements ft^, ^ are independent of i so that, for each 
s, one has d^ copies of the same matrix to diagonalize. 
Once the diagonalizations are performed, one can write 



S d« 2s+l 

^-E EE^i^- 

s=[S] i=l a=l 



(A5) 



where eigenvalues £"1*' of iJ^*' are independent of i and 
where the corresponding eigenvector \s;a)i is given by 



(s) 

with coefficients aa:m G C being independent of i. 
The partition function is then evaluated as 

S rf" 2S-I-1 

= EEE 



= [S] i=l a=l 



E 

s=[S] 



E 

Q = l 



-/3-Bi' 



E '^^^^^^ 

s=[S] 



(A6) 

(A7) 
(A8) 

(A9) 



where Z^^'' — Tre '^^rct is the partition function associ- 
ated to any of the H'f^ (here we choose a reference index 
i — ref). 

In the same vein, one can write the entropy as 



s 



£AB = -Tr[plog2p]= d^£ 

s=[S] 



N p{s) 
AB' 



(AlO) 
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where p 
defined £ 



AB 



Tr 



is the density matrix and where we 



P^cf l0g2 Prof 



and 



Computation of Sa and £b 

We split the system in two subsystems A and B con- 
taining L = tN and N — L = (1 — t)N spins respec- 
tively, with the aim to compute the reduced density ma- 
trices PA = Trsp and pB — Tr^ip. The partial traces 
can be performed most easily by applying the technique 
introduced above, for each subsystem. Therefore, one de- 
composes each spin sector ,s into spin subsector si and S2 
(indices 1 and 2 refer to subsystem A and B respectively) . 
The Hamiltonian then reads 



S Si S2 



^=E E' E'E E4^ 

s=[S] si = [Si] S2 = [S2] il = l i2 = l 



(si,S2;s) 



(All) 



S Sl S2 '^ = 1 "^2 s s 

= E E' E'E E E E 

s=[S] si = [Si] S2 = [5'2] il = l i2 = l m=-s m' = -s 
'^m',™ l*l'S2;S,™')«i,»2 ii,i2(si:S2;S,"l|, (Al2) 

where |si, S2; s, m)i^,i^ denotes an eigenstate of operators 
Sf, S2, and Sz with eigenvalues si(si -I- 1), 52(^2 + 1), 
s(s + 1) and m respectively. Index ii {12) labels the d^_^ 
{d^2^^) degenerate subspaces of spin si (52) that can be 
built from L {N — L) spins 1/2. For a given s, primed 



sums are restricted to values of Si and S2 that can add 
to a total spin s. That is to say, one must fulfill the in- 
equalities |si ~ S2 1 ^ s ^ Sl + S2j so one can alternatively 
write : 

S Sl S2 Sl S2 S1+S2 

E E' E' - E E E • (A13) 

S = [S] Si = [Sl] S2 = [S2] Si = [Sl] S2 = [S2] S=|si-S2| 

We have denoted 5 = N/2, Si = i/2 and S2 ^ (N - 
L)/2 the maximum spins of the whole system and of each 
subsystems. As before, minimum spins are denoted with 
square brackets. Note that the degeneracy of the 
spin-s sector can be recovered from 



Sl , S2 

E E'<^: 

Sl = [Sl] S2 = IS2] 



JV-L _ lAT 



(AM) 



The matrix elements h^J, ^ are the same as in Eq. ( A4) , 
and do not depend on si, S2, ii or Z2. Thus, all Hamil- 
(si,s2,s) Y^^-^^ ^Yie same eigenvalues Ea ^ [which 



tonians 7?, 

are the same as in Eq. (A5)], with the corresponding 
eigenvectors 

s 

|si,S2;s;a)»i,i2 = E ^a^m\si,S2;s,m),^^,^. (A15) 

Once again, coefficients independent of si, 

S2, *i and 12, and are the same as in Eq. (A6). 
Then, the density matrix p — reads 



^ Sl S2 S1+S2 "si 2s+l 

P=^E E E EE E ^"^^"^ I*i''^2;s;a)ii,i2 ^l,^2(sl,S2;s;a|, (A16) 

Si = [Sl] S2 = [S2] S = |si-S2| *1 = 1 *2 = 1 0=1 

= 1 E E E' E E E E E e-^^^^'a(f:a[:^L'ki,-2;s,m'),,, ,,,,(^1,^2;^,™!. (AIT) 

si = [5'i] S2 = [S2] s=|si-S2| ii=l 42 = 1 "=1 m=-s m' = -s 



r 



where a* denotes the complex conjugate of a. 

One of the final steps is to decompose each basis state 
|si, S2; s, m')ij_i2 on the tensor-product state basis using 
the Clebsch-Gordan coefficients 



\si,S2;s,m 

Sl S; 

E 



(A18) 



S2 



nil ,m2 ;m 



/ y ^si,S2;s 

mi= — Sl JTl2 = — S2 



|si,TOl)ii «) |s2,rn2)» 



with obvious notations. The only couples {1111,7712) that 
contribute to the above sum are such that mi +m2 — m, 
since the Clebsch-Gordan coefficients vanish otherwise. 



Introducing the shorthand notation 



Sl S2 S1+S2 31 S2 2S+1 S S Sl S2 Sl S2 

E-EEEEEEEEEEEE- (ai9) 

all si = [Si] S2 = [S2] s=|si— S2I 41 = 1 42 = 1 1=1 m=— s m' = — s mi= — si m2=— S2 t?i'i= — Si m2 = — S2 
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one then gets 

P= ^Ee"''''°'«"f-«il'^™ ^ |s2,m^)., .,(s2,m2|. (A20) 

all 

It is now possible to perform a partial trace. Let us focus on computing pA = Trsp- This partial trace will enforce 

m'2 = m2 in X^aii- Furthermore, the index 12 disappears from the quantity to be summed over, so that '^i^Li simply 
yields a factor d^~^ ■ At the end of the day, one finds 

^ Si "si Si Si 

= ^ 51 Z ^^m'il™il«i'"^i)n with 

si = [Si] ii = l T?ii= — si ni'i= — si 

S2 S1+S2 max^ji /2s+l \ 

^ _ rl^-L \ " ( (s)*„(«) 1 ^mi,m-mi;m ^mi,m-™i;m+™'i-mi /Aon 

'^'4„'^^„ii - 2^ "s2 I '^a;m"a;m+m'i-mi I '-^si,S2;s L'si,S2;s ; l^/ij 

S2 = [S2] s=|si — S2I m=minm \ a=l / 



with 

mium = max(— s, —s + mi — m'^, — S2 + m-i), 
max,„ = min(s, s + toi — to^, S2 + fni). (A22) 

These limits come from the fact that all values 
of m are not allowed in between —s and s, but 
one must also ensure that —S2 ^ m — mi ^ S2 and 
—s ^ TO + m'l — TOi ^ s. Let us remark that, here, we 
have got rid of the sum over m2 and kept the sum over to. 
One could have of course done the opposite, resulting in 
a somewhat different implementation. Expression (A21) 
of Pa allows one to compute £a- The partial entropy £b 
can be computed similarly. Once these are known, the 
mutual information T is simply obtained from its defini- 
tion (9). 

It should be noted that the above calculation requires 
the numerical computation of many Clebsch-Gordan co- 
efficients. As there are too many of them, these cannot 
be stored and must be computed when they are needed. 
This was achieved by implementing the algorithm devel- 
oped by Schulten and Gordon in Ref. 48, which is fast 
and stable. 

As a final comment, let us give an estimate of the 
computational complexity of our algorithm. Looking at 
Eq. (A21), one sees that eight sums are involved. How- 
ever, the sum over ii will only yield a multiplicity, which 



leaves seven sums. Furthermore, the values of the inner 
sum over a, namely 

E (A23) 

Q = l 

can be computed once (for a given s), and stored in a 
matrix. This matrix is nothing but the density matrix 
in the spin-s sector (up to a factor of 1/^). The 
computation cost of this sum is thus negligible, and one 
is left with a sum over six variables, namely si, mi, 
to'i, S2, s and to, which all take a number of values 
scaling linearly with N (at most). So we expect that our 
algorithm roughly scales as (which we observed in 
practice) . 
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